; docformat = 'rst'
;
; NAME:
;   MR_MAPDIST
;
;+
;
; :Purpose:
;   Calculate 2D array of distances to grid cells from a particular point
;
; :Categories:
;   Model utilities
;
; :Params:
;   LonC : Longitude of point to calculate distances from
;   LatC : Latitude of point to calculate distances from
;   Lon : Longitude array (1D)
;   Lat : Latitudes array (1D)
;   
; :Keywords:
;   NoWrap : Set to TRUE to prevent wrapping around edge of map (only on longitude direction)
;
; :Examples:
;   Calculate array of distances from 10N 20E::
;     IDL> distances=mr_MapDist(10., 20., lon, lat)
;
; :History:
;   Written by: Matt Rigby, MIT, May 31st 2011
;
;-
Function mr_mapdist, lonc, latc, lon, lat, nowrap=nowrap

  lon=double(lon)
  lat=double(lat)
  lonC=double(lonC)
  latC=double(latC)

  lonSize=n_elements(lon)
  latSize=n_elements(lat)

  mapdist=dblarr(lonSize, latsize)

  For LonI=0, LonSize-1 do mapdist[LonI, *]=sqrt( (lon[LonI] - lonC)^2. + (lat - latC)^2. )
  ;wrap around distances
  if keyword_set(nowrap) eq 0 then begin
    if lonC gt 180. then begin
      wh=where(lon le lonc-180.)
      for lonI=wh[0], wh[n_elements(wh)-1] do mapdist[lonI, *]=sqrt( (lon[lonI] + (360. - lonC) )^2. + (lat - latC)^2.)
    Endif
    if lonC lt 180. then begin
      wh=where(lon ge lonc+180.)
      for lonI=wh[0], wh[n_elements(wh)-1] do mapdist[lonI, *]=sqrt( (lonC + (360. - lon[lonI]))^2. + (lat - latC)^2.)
    endif
  endif
  
  return, mapdist

End